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Abstract 

In this paper a new approach is derived in the context of shape theory. The im¬ 
plemented methodology is motivated in an open problem proposed in |Goodall and| 
Mardia (1993) about the construction of certain shape density involving Euler 


hypergeometric functions of matrix arguments. 

The associated distribution is obtained by establishing a connection between the 
required shape invariants and a known result on squared canonical correlations avail¬ 
able since 1963; as usual in statistical shape theory and the addressed result, the 
densities are expressed in terms of infinite series of zonal polynomials which involves 
considerable difficulties in inference. Then the work proceeds to solve analytically 
the problem of computation by using the Eulerian matrix relation of two matrix 
argument for deriving the corresponding polynomial distribution in certain para¬ 
metric space which allows to perform exact inference based on exact distributions 
characterized for polynomials of very low degree. A methodology for comparing 
correlation shape structure is proposed and applied in handwritten differentiation. 


1 Introduction 

Invariant statistics has emerged as one of the most powerful tools in statistics to study ma¬ 
trix variate problems, in particular the statistical shape analysis has been rolled a number 
of applications in almost all the natural and exact sciences; their applications in machine 
vision and images analysis are notable. In general, shape theory deals with problems of 
differentiation of “objects” which can be summarized by anatomical or mathematical marks 
called landmarks, under the assumption that those points can be placed (in presence of ran¬ 
domness) in any similar objects and absorbs all the geometrical information of the “shape” 
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of the object. The invariant theory enters in the analysis, when the comparisons, match¬ 
ing, and/or inference are performed in special quotient spaces rather than the usual full 
noised Euclidean space. Then according to the symmetries involved in the study defined 
and supervised by the expert in the experiment, a number of different shapes can be defined. 
Three classical shapes are well studied in literature, the Euclidean or similarity shape, the 
affine or configuration, and the projective shape; then the associated shape is all geometri¬ 
cal information which remains under removing location, rotation, scaling, reflection, uniform 


shear and/or projection, etc., see for example, Dryden and Mardia (1998) and the references 
therein. 

Under deterministic assumptions, shape theory counts a number of approaches, a source 


of their main aspects are given for example in Kendall et al. (1999), their references; as 


a consequence of this geometrical understanding of shape, the procrustes theory, including 
projective shape analysis and its relation to affine shape analysis, have been studied profusely 


in applications; some works in these directions are due tc 

Mardia et al. 

(200, 

3), Patrangenaru 

and Mardia 

2003) 

Mai et c 

il. (2011), 

Domokos and Kato (2010), 

Ecaberth and Thiran 

(20041, 

Glasbey and Mardia ( 

2001) 

Goodall (1991J, Groisser and Tagare ( 

2009), Lindeberg 

and Carding 

1997) 

Mardia et al. ( 

1996) 

Mardia and Patrangenaru (2 

005), 

Mokhtarian and 

Abbasi 

(2002 

, Horgan et al. 

1992 

), Kent et al. (2004), 

Lin and Fan^ 

(2007), among many 


others 

Now, the statistical theory of shape deals with landmarks in presence of random (see for 
example Dryden and Mardia (1998), Small (1996) and the references therein). A branch 
of the statistical shape analysis was connected by Goodall and Mardia ( 1993) (corrected 


by Dfaz-Garcfa et cd\ (2003) and Caro-Lopera et al. (20101)) with the very old theory of 


matrix variate distribution based on A.T. James zonal polynomials (see for example James 


(19641 and the cited papers). At that time the shape distributions derived in Goodall and 


Mardia (1993) and most of the classical results in matrix variate distributions collected 


in the excellent book of Muirhead (1982) could not computed and force to use asymptotic 


distributions. Only recent the zonal polynomial were computable in an efficient way by|Koev 


and Edelman (20061, and a hope for computations of infinite series of zonal polynomials was 


a possible task under certain restriction and subject to the problem of automatic detection 
of convergence of the series. 

Parallel to the problem of computation, the classical shape theory of |Goodall and Mardia| 
(1993) based on the Gaussian distributions, opened a field of new interesting theoretical 


and applied problems, the so called generalized shape theory, which models shape under 


non Gaussian models, see for example, Caro-Lopera et al. (2009b 2010); Caro-Lopera et al. 


(2014c 

), Dfaz-Garcfa & Caro-Lopera 

(2012a|b 

20141; 

Dfaz-Garcfa and Caro-Lopera 

(2016) 

and 

Caro-Lopera et al. 

(2013a 

), among many others. Moreover the generalized shape anal- 


ysis can be studied in the more general setting of real normed division algebras, a unified 
approach valid for real, complex, cuaternion and octonion landmark data matrices Dfaz- 
Garcfa & Caro-Lopera (2010). 


As is can be check in the above references the Euclidean shape based on similarity trans¬ 
formations, in all theirs forms, and the affine configuration have been studied extensively, 
involving inference, families of distributions and so on; however, an open problem seems to 


be unsolved from the times of Goodall and Mardia (1993). 


Analysis of shape under Gaussian models were proposed by Goodall and Mardia (1993) 


following a natural induction way of deriving new shape distribution via Laplace-transform in 
order to obtain hypergeometric series of matrix arguments by known integration over cones 
and manifolds involving zonal polynomials. This methodology was though for obtaining 
successive distributions which filters new more valuable information of the original landmark 
data. 

The underlying sequence of new shape distributions appeared as follows: start with an 
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original Gaussian landmark data, then after removing traslation and rotation the so called 
size and shape distribution is obtained, this function has the form of a qF-\ (•) hypergeo¬ 
metric function Bessel type (see Muirhead] (1982)), a known series of zonal polynomials. 
According to the inductive process of constructing hypergeometric series, the next appro¬ 
priate shape should be a confluent type i-Fi(-) by removing the scaling information of the 
size and shape, then the shape distribution is gotten. But as it can be seen in |Goodall and| 
(1993), this shape has not the required form and no properties were available at 


Mardia 


that time to explore the asymptotics or computation of the distribution (later Diaz-Garcfa & 
Caro-Lopera (2012a|b 2014); Diaz-Garcfa and Caro-Lopera (20161 generalized the result to 


every elliptical model and transformations and performed inference with the exact densities). 
So, in order to obtain confluent shape distribution, Goodall and Mardia | (1993) defined an 
affine transformation which lead directly from the Helmertized matrix (the original Gaus¬ 
sian landmark matrix invariant under traslation) to the required confluent distribution (later 
Caro-Lopera et al. ( 2009b| 2010); Caro-Lopera et al. (2014c) generalized it to the class of 


elliptical distribution and performed inference with the exact distribution, some of them 


simple polynomials of low degree). Finally, Goodall and Mardia (1993) conjectured that 


a new class of shape distributions could be constructed by using the induction way and 
claimed for an hypergeometric type shape distribution 2 Fj(-). 

In this paper we study the addressed open problem by defining a new approach in the 
context of statistical shape theory. Section [2] focus in the definition of the new shape in¬ 
variant. As very usual in shape theory, the resulting density is a series of zonal polynomials 
which presents strong problems for inference; it motivated Section [3] which solves the prob¬ 
lem of inference by considering in some simple parametric space an equivalent polynomial 
distribution of very low degree in the latent roots of the matrix arguments of the zonal 
polynomials, then the exact formulae available for these polynomials in two dimensions, a 
number of applications and methodologies for shape comparisons are given in Section [4] 


2 Eulerian Shape 


The aim of this approach consists of proposing an alternative way to compare populations 
in the context of the shape theory. It is easy to check in the literature shape that Wilks’ 
theorem is the most common resource for testing means of two populations (a number of 
other tests can be implemented, see Dfaz-Garcfa and Caro-Lopera (2008)), |Caro-Lopera 


et al. (2010 ?). However, the application of the classical chi-squared distribution of that 


theorem demands a number of conditions on regularity which are difficult to keep in most 
of the available landmark data. 

The strength of relation between two variates has a similar power of explanation in shape 
theory than the equality of means and both complement each other for completing a shape 
analysis. 

Consider two NxK matrices X** and Y**, which summarize two figures comprised in N 


landmarks in I\ dimensions, such that X** ~ A/’(/x- ? 


L JV l 




X* 


.) and Y** ~ A 




3>y**)j where <I? X »* and 3 ?y** are positive definite K x K matrices. 

Now, as usual in shape theory we remove some meaningless geometrical information as 
the translation of both figures, by using contrast matrices, for example sub-Helmert arrays. 
In this case we use a (N — 1) x N matrix, L, which has orthonormal rows to 1 = (1, ..., 1)'; 
to construct the so called Helmertized N — 1 x K matrices, say, X* and Y*, respectively, 
via LX** = X* and LY** = Y*. 

Under these transformations, the randomness of the transformed landmarks follow the 


laws: X* ~ N(n 


Y* > E ® 5jh), where n = N — 1 x K . The means 
are changed to ^t x „ = L/i x »» and = L/i, Y «.», while the new covariance matrices are 


X* ) *-n ' 


>s 22 ), Y*~AA(m 
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conveniently denoted by Su and S 22 , where Sn = L3 >y**L' and S 22 = L3>x**L'. 

We storage the above Helmertized matrices into an N x 2K partition array W = [Y* |X*], 
which implies that W ~ A/”([/x Y * |Mx*L <8> S), where 

■y _ ( Su Sl2 

_ ^ S 21 S 22 

Finally, set Z = W — [/n Y . |/x x *]> thus Z = [Y|X] ~ Af( 0, I n ® X), where we have defined 
Y = Y* — /Lt Y * and Y = Y* — /li y » • 

Now, we are interested in the measure of the relation between the two variates X and Y, 
where X 12 = X' 21 is the covariance matrix between the components of X and the components 
of Y. In order to establish that comparison we need to define the following supporting 
equivalent class of shapes. 

Definition 2.1. Given two Helmertized n x s matrices X and Y, it is said that they have 
the same eulerian shape, if 

1 . hx = [x; 0 ]', where Xi is a non-singular s x s matrix and H £ 0(r), 

2. HY = T = [U'|V']', where U is s x s, V is (r — s) x s and H is defined according to 
(1), and, 

3. U and V satisfies |(U'U)(U'U + V'V) -1 - I s | = 0. 

In this definition note that the matrix X is of rank s, then there exists a rotation matrix 
H £ 0(n) in order that HX = [X) 0]', where Xi is a non-singular s x s matrix, this matrix 
absorbs the eulerian meaningful shape of X, and then the matrix Y is rotated with that 
preceding orthogonal matrix, and partitioned according to HY = T = [U'|V']', which must 
satisfy |(U'U)(U'U + V'V) -1 — I s | =0. Observe also that U and V depend on X. 

The name of this class of shape obeys the type 2 Fj(-) of the hypergeometric series involved 
in the main theorem of the present work. 

We can put together the above transformations and definitions into the following steps 
which obtains the eulerian coordinates (U, V) of a given pair of original figures (X**, Y**): 


L ( Y** X**) = (Y* X* ) = W = ( Y + My , X + /x x » ) => 

subject to |(U'U)(U'U + V'V) -1 - I s | = 0. 

Geometrically, Helmertized matrices means that Y* and X* are invariant under trans¬ 
lation, and then with the remaining transformations, the resulting eulerian shape (U, V) 
of the pair (X*, Y*) is all geometrical information remained after filtering out the rotation 
and permutations of the columns from the Helmertized pair (X*, Y*). 

This last property can be understood as the eulerian shape is invariant under the label 
of the axis of coordinates. For example, under this equivalent class, figures in an Euclidian 
three dimension space labeled with points with coordinates (x, y, z) are equivalent to figures 
with permuted the axis of coordinates, as (y, z, x), z, x, y , (y, x, z), (z, y, x) and (x, z , y). This 
unusual invariance does not hold in the classical shapes via euclidian, affine or projective 
transformations and it is the most important novelty of the proposed approach, because does 
not require that the landmarks of the two objects strictly match in labels of the coordinates 
of axis. 

The few works published in this area belong to procrustes methods based on unlabelled 
points where the minimizations of the distances between the figures are carried out over the 
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similarity transformations and permutations of the landmarks labels, they include a number 
of applications as megalithic standing stones ( Cornwall) and human vision, among others, 
see section 12.9 of Dryden and Mardia (1998). As usual, the problems of inference and 


the exact density are one of the obstacles to solve in the existing approaches. Finally, this 
invariant statistics of shape analysis enters too into the so called landmark free approaches, 


a wider area of study, see for example section 12.10 of Dryden and Mardia (1998). 


Clearly, with definition [2d] we can implement classical procrustes analysis, as it was done 
for Eucliand and affine transformation, but we are interesting in the statistical analysis of 
the emergent shape. 

In presence of randomness in the location of landmarks in the two Helemertized figures 
X and Y (linear transformations of the original figures X** and Y**), the so called squared 
canonical correlation coefficients allow the statistical matching of the objects via 

= 0 , ( 1 ) 


(Y'Y)- 1 Y'X(X'X)- 1 X'Y - r 2 I K \ 


which is reduced to 


|(U'U)(U'U- 


V'V )- 1 — r 2 I 


K 


= 0 , 


( 2 ) 


by the items (1) and (2) of Definition 2.1 

As we proved, X and Y are the n = N — 1 x K helmertized matrices, such that Y ~ 
7V(0,I n <g> Eh), X ~ A/"(0, I n (g) S 22 ) and we are interested in the measure of the relation 
between X and Y, where E 12 = S 21 is the covariance matrix between the components of 
X and the components of Y. Thus, Z = [Y|X] ~ A/”(0, I n 0 E), where 

E11 E12 

X21 X22 


E = 


The invariant statistics involved in the shape theory context arrives here when we study 
the distribution of the so called squared canonical correlation coefficients under the Eulerian 
transformation above described, which obtains the squared canonical form of the covariance 
matrix as: 

/IPO 
S= p I 0 
\ 0 0 I 

where P is the K x K diagonal matrix with the squared canonical correlations pi,..., pi< 


arranged in decreasing order down the diagonal (see Constantine (19631, Muirhead (1982), 
p. 539). 


So, from theorem 11.3.2 of Muirhead (19821, p. 557 (see also Constantine (1963)), we 
have: 

Let A = Z'Z = ^ x'Y X'X ) = ( A^ A^ ) ^ ave the W 2 x(n, E) distribution, 

where n > 2K. Then the joint probability density function of r\,... ,r 2 K , the latent roots of 
A- 11 1 Ai 2 A 2 2 1 A 2 i, is 


K 


f(rl 


: 1 K 


) = n 


(K) 


In, In-ll<- P 2 ,R 2 
2 2 2 


t k 2 /2 


r K (in) 


T \{\K) T k [\{n-I <)] 

K 


X n [(>-,T 1,2 (l - r?)<»-“'-■>/*] - r% 


i<j 


for 1 > r\ > - - - > ri > 0, where p\,...,p 2 K are the latent roots of E 11 1 Ei 2 S 22 1 E 
P 2 = diag(p 2 , ...,P 2 k ) and R 2 = diag(r^,..., r 2 K ). 


21 ) 
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3 The main result 


For decades, the above type of expressions could not computed and they forced to asymptotic 


studies (see section 11.3.5 ol 

Muirhead 

(1982 

)), and it seems that with the recent algorithms 

of Koev and Edelman 

(2006| the exact distribution can be used for applications. However, 

recently 

Caro-Lopera et al. 

(2014c 

1, proved that those algoritms do not converge when the 


inference is performed in certain shape invariants as the affine case even applied to small 
size landmark data. 

Computation in matrix variate theory via zonal polynomials (see Muirhead (1982)) and 


invariant polynomials of several matrix arguments extending the preceding ones (Davis 


(1979), Davis (1980)) has been the biggest problem in such area, there are hundreds of 


papers deriving numbers of theoretical properties involving such polynomials or asymptotic 
distributions, etc. but a very few, considering the problem of inference with the exact 
densities. For decades the zonal polynomials were impossible to compute until the works of 
Koev and Edelman (2006), meanwhile the computation of invariant polynomials of Davis 


remains unsolved since 1979. However, infinite series of zonal polynomials has considerable 
problems and so the performance of inference associated to the large number of results. 

The problem for the computation of one hypergeometric series is addressed by Koev] 
and Edelman (2006) in page 845 : “Several problems remain open, among them automatic 
detection of convergence”. This difficulty is amplified when we have to consider the opti¬ 
mization of a large product of such hypergeometric series, and the addressed algorithms do 
not allow large truncations to obtain certain accuracy and the convergence is completely 
doubtful. 

So, in this scenario we have two ways to follow, a computational study which involves 
also a machine problem or trying an analytical solution. 

Motivated by a similar situation in affine shape studied by Caro-Lopera et al. (2014c), 
this section attends an analytical solution to the problem of Eulerian shape inference. As 
before, the solution resides in the fact that shape densities are usually infinite series of zonal 
polynomials, but a few of them can be studied in the setting of the so called polynomial 
shape densities, which exploits some characteristics of the parameters of the distributions 
indexed by the number of landmarks and dimensions, and by using matrix generalized ma¬ 
trix relations, the density can be turned into a polynomial. For example, the affine shape 
theory is characterized by a confluent hypergeometric function type which can be turned 
into a polynomial density by some special domain and using the so called generalized matrix 
Kummer relation of Diaz-Garcia and Caro-Lopera (2015), a generalization of the well known 
Kummer relation by Herz (1955) with gaussian kernel. Continuing in this way, we expect 
that the next type of hypergeometric shape distribution, the Eulerian of this paper could 
be transformed into a polynomial by using some special parametric space and the similar 
Euler matrix relations. However, the density derived in the previous section is a modified 
Eulerian hypergeometric type and involves two matrix arguments instead the expected one, 
which facilitates the automatic application of the Euler matrix relations. Given the con¬ 
nection between the Eulerian shape density and the squared canonical correlation result 
of Constantine (1963), it seems that non existence of that relations explains why the last 
result was set as an infinite series of products of zonal polynomials, instead as a polynomial 
of low degree, as it happens with a number of result based on Gaussian distributions and 
Kummer relations. Then an improvement of the shape density will be applied in the squared 
canonical correlation theory too. 

Specifically, recall that the hypergeometric function of two m x m matrices X and Y as 
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arguments are given by 


D i* m >(a 1 ,...,a p ;& 1 ,...,&,;X, Y) = £ £ 


k—0 k 


(ai) K ---(a p ) K C K (X)C K (Y) 

(for )*•■■(&,)* fc!C K (I m ) ■ 


( 3 ) 


The problem is reduced to obtain certain matrix relation which turns certain series of the 
type 2 f{ m ' > (a, b; c; X, Y) into an expression of the type /(c — a, 6; c; X, Y) in such way that 
c— a are negative integers or negative half-integers, which will turn /(•) into a polynomial 
of low degree. 


The required computable is obtained from the following matrix relation due to Caro- 


Lopera (2014): 


2^(0!, 025 6 i; X,Y ) = EtiE 

k =0 K ' k 


1 [b\ — a 1 ) K (a 2 ) K 


(bi)h 


( 4 ) 


>0(r, 


|I m - XHYH' \~ a2 C K (—XHYH'(I TO - XHYH')“ 1 )(dH); 


if bi — a\ is a negative integer, say —q, the hypergeometric function is a polynomial of degree 
mq. 

And given that n = N— 1(> 2K) is the number of helmertized landmarks in K dimension, 
then f — § = (K — N + l)/2 is a negative integer or a half negative integer, so we have 
arrive to main result of this work: 

Theorem 3.1. Let A = Z'Z = ^ x'X ) = ( A U A^ ) /ui,Je the w 2 if(n,S) 

distribution, where n > 2K. Then the joint probability density function of rf ,..., r\, the 
latent roots of Aj" 1 1 A 12 A^ 2 1 A 2 i, is 


K 


K 


f(rl 


,rr k i\ K ) [b(n- k)\ n^ 1 rf )B/a nw ^ 


if 


n [Mr ,/2 d - 4) < -“- i>/2 ] x s x 

■ ' 1=0 A 

2tt/ I —n / 2 


( A _ (A) 

V 2 2 > \\2 I) 


i=1 


x / I/f — P HR H 

io(jr) 


(K) 

1=0 ” A V 2 / A 

Ca (-P 2 HR 2 H'(I^ - P 2 HR 2 H') _1 ) (dH), 


for 1 > rj > • • • > rf > 0, where p\,..., p 2 K are the latent roots of ^11 ’^12^22 ^21, 
P 2 = diag(p 2 ,..., p^-) and R 2 = diag(rf ..., r^). // f — ™ is a negative integer, say —q, 
the function is a polynomial of degree Kq in K variables. 


Note that the parametric condition of Eulerian polynomial shape, i.e. f (n > 2K) 
just demands that the objects in AT-even(odd) dimensions are summarized by a number of 
TV-odd(even) anatomical landmarks. 

Nevertheless the most of the literature of statistical shape analysis is valid for any dimen¬ 
sion K , the application of the theory seems to be very restricted to anatomical landmarks 
placed on two and three dimensional objects, in fact the classical text on shape theory, such 
Dryden and Mardia (19981 only uses 2-D landmark data, and 3-D applications are very 


as 


rarely given the complexity of the distributions and the registration of invariant critical 
points. So, when we consider K = 2, all the typical landmark data presented for example 


Dryden and Mardia (1998) can be studied and the integral over the orthogonal group 
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is a feasible and programming task using the exact formulae for Jack polynomials (which 
includes real zonal polynomials) due to Caro-Lopera et al (2007). 

Next we give some summary of Jack polynomials and exact two dimensional formulae, 
which particularized to zonal polynomials, will be needed in some applications. 

Let us characterize the Jack symmetric function J^ a \yi, ■ ■ ■ ,y m ) of parameter a , see 


Sawyer (1997). Consider a decreasing sequence of nonnegative integers k = (k\,k 2 ,. ■ 


with only finitely many nonzero terms is said to be a partition of k = Let n and 

A = (Zi, Z 2 , - - -) be two partitions of k. We write A < n if 1 l i < 1 k i for each t. 
The conjugate of k is n' = ...) where k\ = cardjj : kj > i}. The length of n is 

l(k) = max{i : fc,; ^ 0} = k[. If 1(k) < in, one often writes k = (fci, k 2 , ■ ■ ■, k m ). The 
partition (1,..., 1) of length m will be denoted by l m . 

And recall that the monomial symmetric function M K (-) indexed by a partition n can be 
regarded as a function of an arbitrary number of variables such that all but a finite number 
are equal to 0: if yi = 0 for i > m > 1 (k) then M K (yi ,..., y m ) = y ° 1 ■■■y% n , where 
the sum is over all distinct permutations {<ti, ..., a m j of {ki,..., k m }, and if 1(k) > m 
then M k (y-\,, y m ) = 0. A symmetric function / is a linear combination of monomial 
symmetric functions. If / is a symmetric function then f(y ±,..., y m , 0) = f{yi, ..., y m )- 
For each m > 1, f(yi ,..., y m ) is a symmetric polynomial in m variables. 

Thus the Jack symmetric function ■ ■ ■ ,y m ) with a parameter a, satisfy the fol¬ 

lowing conditions: 


J[ a \yi,---,ym) = '52j K ,\M x (y 1 ,...,y m ), 

m — i + 1 N 


\<K 


t i V a / fc, 

^ 2 d 2 ji a \yi,... ,y m ) 1 dJ^ a \yi, ■ ■ ■ ,y m ) 

2^ Vi Ft,, 2 + n, 2^ y i 2^ 


(5) 

( 6 ) 


dy\ 


a 


s y* - 


dyi 


V ki(ki - 1 + -(m- i))Ji a \y i,■ ■ ■ ,2/m)- 

/ rv 


(7) 


*=l 


Here the constants j Kj a do not dependent on y[s but on k and A, and (a) n = nr=i( a + 
i— 1). Note that if m < 1 (k) then ji a \yi, ■ ■ ■ ,y m ) = 0. The conditions include the 
case a = 0 and then ji°\yi,.. ■ ,y m ) = e K > ~ * + l ) ki , where e K (yi,... ,y m ) = 

n!=i e ki (yi, ■ ■ ■ 1 2/m) are the elementary symmetric functions indexed by partitions k, if m > 
1 (k ) then e r {yi,. .,y m ) = E,j<i 2 <-<v 2/*i' "fc. and if m < 1 (k) then e r {yi, ■. .,2/m) = 0, 
Sawyer] (19 971). _ 


see 


Now, from 


Koev and Edelman 


(2006), the Jack functions jjf‘\ Y) = ji a \yi,... ,y m ), 


with yi,, y m being the eigenvalues of the matrix Y, can be normalized in such a way that 

Y) = (tr(Y)) fe , 


where C“(Y) denotes the Jack polynomials. They are related to the Jack functions by 

C a A Y) = — JT(Y), 

J K 


( 8 ) 


3k = ]J 


(bi)eK 


where 

























and h*(i,j) = kj — i + a{kt — j + 1) and h* K (i,j) = kj 
lower hook lengths at (i,j) € k, respectively. 

Then by applying ([8]), we can write ([7]) as 




2 .2 1 9Ci a) ( Y) 

fg y i-Vj % 


1 + a(ki — j) are the upper and 


5>(fe-l+-(m-0)<?£*)(Y). 

^' rv 


(9) 


i= 1 


Now, when m = 2 in ®>l |Caro-Lopera et a7~~| ( |2007[ ) found the following formulae for Jack 
Polynomials of the Second Order 


0 (fci,fc 2 ) 

r< (“) 
°(fel,fc2) 


00 

(la) 


(2/il/2) (fcl+fc2)/2 AiT 


/ P P 1.1. (i/1 +l/2) 2 \ 
V 2’2 a’2’ 42 / 12/2 / 


(piP2) (fcl+fc2 ~ 1)/2 , F /1 1 + P 1_P 3 (pi+2/ 2 ) 2 \ 

2(2/i + 2/2)- 1 2 U 2 ’2 2 ’ 2 ’ 4 /A 2/2 y 1 ’ 


( 10 ) 


with p being either even or odd. For distinguishing the case under consideration, odd or even, we 
will use the upper indices o or e with j4i and A 2 . Then the corresponding solutions are the following 
Even case. If p = fci — fe = 2n, n = 0,1,2,... then 


A 


e 

1 


(-ir n( i+2 *) 


i=0 



and A\ = 0. 


Odd case. If p = ki — *» 2n + 1, n = 0,1,2,... then 


A° = 0 and A% = (2 n + 1)A1. 


Three particular cases are of interest in the literature: the quaternionic case (a = 1/2), the 
complex zonal polynomials (a = 1) and the real zonal polynomials (a = 2), these results are 
summarized in the following table: 
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a 

P 

a 

b 

c 

Ai 

A.2 

1 

2 

even 

—n 

n + 2 

1 

2 

(—l) n 3 

(2n + l)(2n + 3) 

0 


odd 

n + 3 

—n 

3 

2 

0 

(—l) n 3 
(2 n + 3) 

1 

even 

—n 

n + 1 
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2 

(-1)" 

(2n + 1) 

0 


odd 

n + 2 

—n 

3 

2 

0 

(-1)" 

2 

even 

—n 

n+ 1/2 
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2 

(—l)"(2n)! 
2 2rl (n!) 2 

0 


odd 

n + 3/2 

—n 

3 

2 

0 

(—l)"(2n + l)! 
2 2n (n!) 2 



4 Applications 


In this section we describe the general procedure for applications of the Eulerian shape distributions 
and compare the computational performance when the inference of a certain appropriate landmark 
data for this data is studied under the modified classical result of Constantine {19631 via the 


algorithms of Koev and Edelman (2006) and the exact polynomial distribution here derived. 

For the first aspect, there are a number of situations to consider, for example, the Eulerian 
shape can be used to establish correlation structure of two populations, or determine the correlation 
structure of a population and a given template, also, the correlation structure of a population and 
a given object probably belonging to another population, and finally, we can use the polynomial 
distributions for exploring the problem of landmark discrimination. Next we explain the associated 
procedure of each one. 


4.1 Correlation structure of two populations 

Consider a sample of m figures of a population X and m figures of a population Y; both samples 
satisfying the shape invariance assumed for the Eulerian shape and the fact that (K — N + l)/2, 
(N — 1 > 2I \) is a negative integer, i.e. the objects in A'-even(odd) dimensions are summarized by 
a number of IV-odd(even) anatomical landmarks. For each pair of figures Z = [Y|X] compute the 


corresponding rf j. 
the function fj(rf, 


= 1 ,.. 


’ r KJ, J 

,ri) of Theorem 3.1 


m the latent roots of A^A 


A„ 1 A 2 i, in order to obtain 


12^-22 


10 














































Then maximize the likelihood function 


3 = 1 

respect to the population squared canonical correlation coefficients pi,..., p 2 K . 

These estimates will provide the strength of relationship between the two populations in the 
sense of the Eulerian shape invariance. 

As we shall see in the experiments various methods can be implemented in order to obtain more 
information about the shape equality of both populations without using the classical statistics for 
means, instead we can use the exact distribution of certain template vs one of the population and 
obtain the probability than an extreme correlation structure of the other population influenced by 
the same template, can be occurred. 

4.2 Correlation structure of a population and a given template 

Exactly as before, but in this case, we have m repeated templates X and we want to verify the 
correlation between the population Y and the template X under the Eulerian shape assumptions. 

4.3 Correlation structure of a population and a given figure 

Exactly as before, but in this case, we have m repeated figures X and we want to verify the correla¬ 
tion between the population Y and a given “figure” X, i.e. it can explain, under the requirements 
of the Eulerian invariants, if the last one belongs or not to the population Y. 

4.4 Landmark discrimination 

This is one of the main less studied problems in the shape theory. It concerns the reduction of 
set of landmarks for a given figure without changing the explanation of the object. In other words 
it worries about reduction of dimensionality of a figure and set the shape theory in the domain 
of very low dimensions. The method of the preceding section can be used in order to study the 
landmark discrimination. To get this end, we start with figures Y and X (in any of the three 
possible scenarios explained in the above items) with a large number of landmarks N (not necessary 
anatomical landmarks, they can be mathematical or pseudolandmarks), estimate the corresponding 
p's, then reduce the number of landmarks and repeat the estimation until the p's change drastically, 
then we reach an optimal number of landmarks which summaries properly the object. 

This technique is particular appropriate for shape analysis involving mathematical or pseu- 
domathematical landmarks, which usually are equally-placed in the contour of the object. For 
example, the method is applicable in some extreme cases, such as discrimination of classes of pota¬ 
toes or differentiation of non regular objects with no chance of placing anatomical landmarks objets. 
A profuse study of this landmark discrimination methodology is will be carried out in a subsequent 
work. Another useful application under consideration tries an automatic detection of anatomical 
landmarks using their limit from pseudo critical points. For example, the 46 pseudo-landmarks de¬ 
fined by Dryden and Mardia ( 1998| in the contour of the mouse vertebrae can be used for studying 
the discrimination problem and finding the possible optimal set of six anatomical landmarks usually 
referred in literature. 


4.5 An application in postcode recognition 

We end this section by illustrating the addressed problem of inference with the distribution provided 


by Constantine (19631 and the algorithms for hypergeometric functions by Koev and Edelman 
(2006). Then we compare the above results with the polynomial density here derived. 

We select a known landmark data of literature on postcode recognition which has presented 
problems in inference using another shape invariants as the affine shape, see |Caro-Lopera et al.| 
(2014c). 
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In this case, we have a 30 sample of handwritten digit 3 (see Figure [TJ, with N = 13, ( I\ = 2) 
landmarks placed according an equally spaced two joint octagons as a template (see Figure |2|. We 
are interested in estimating the strength the estimate of correlation structure between the template 
and the handwritten digit 3. The main theorem can be used to obtain the estimates of population 
squared canonical correlations pi and p% by optimizing the corresponding likelihood. In this case 
the distribution is a polynomial of degree — \{K — N + 1) x K = 10 in the two latent roots of the 
matrix argument in the zonal polynomial. 


oo 

T5 

I 

8 30 


Sample of 30 handwritten digit 3 


..-**••* • » « • 

* » • 

* * « * * 
* "*•:'«,* • » * 

« « » ►. . * . *'» » * : * 

* • ***'•*'•* 
* » X 

* * • 

*. ■ »•' * .■'*.•• * 


' *, * ■ * * 
* * * » * * * 

* ***** * 
* *-* » * 
* •*,,.» , * . « * * 
* * * * -* * « '.* 


• '*.*,*•* ' « * * . * 
.**..» * 
* * * • » 


* * * 


* * * 


10 15 


Original Euclidean x-coordinates of the sample 


Figure 1: 30 handwritten digit 3. 


If we truncate the origi nal distributions derived by Constantine (19631 and use the algorithms of 
Koev and Edelman (20061 we obtain the following table, which shows the starting point (r^, r^)for 
the optimization, the population estimate coefficients between the template and the handwritten 
digit 3, the convergence and the number of iteractions of the algorithm. In this case we just truncate 
the infinite series of zonal polynomials at 10th degree in order to compare the performance with the 
polynomial of 10-th degree which provides the solution. The computations were performed with a 
processor Intel(R) Corel(TM)2 Duo CPU, E7400@2.80GHz, and 2,96GB of RAM. 



r 2i 

Fval. 

Conv. 

Iter. 

pi 

P2 

0.1 

0.1 

372,8472 

y 

23 

0.1260 

0.1260 

0.1 

0.2 

372.8472 

y 

23 

0.1260 

0.1260 

0.1 

0.3 

372.8472 

y 

27 

0.1260 

0.1260 

0.1 

0.4 

372.8472 

y 

30 

0.1260 

0.1259 

0.1 

0.5 

372.8472 

y 

44 

0.1260 

0.1260 

0.1 

0.6 

372.8472 

y 

31 

0.1260 

0.1260 

0.1 

0.7 

372.8472 

y 

31 

0.1260 

0.1260 

0.1 

0.8 

372.8472 

y 

32 

0.1260 

0.1260 

0.1 

0.9 

- 

n 

- 

- 

- 



r 2i 

Fval. 

Conv. 

Iter. 

pi 

P2 

0.2 

0.2 

372.8472 

y 

27 

0.1260 

0.1260 

0.2 

0.3 

372.8472 

y 

34 

0.1259 

0.1260 

0.2 

0.4 

372.8472 

y 

32 

0.1259 

0.1259 

0.2 

0.5 

372.8472 

y 

37 

0.1259 

0.1260 

0.2 

0.6 

372.8472 

y 

36 

0.1260 

0.1260 

0.2 

0.7 

- 

n 

- 

- 

- 

0.2 

0.8 

- 

n 

- 

- 

- 

0.2 

0.9 

- 

n 

- 

- 

- 


r ii 

r 2 

' 2i 

Fval. 

Conv. 

Iter. 

Pi 

■s> 

to K 

0.3 

0.3 

372.8472 

y 

29 

0.1260 

0.1260 

0.3 

0.4 

372.8472 

y 

39 

0.1260 

0.1260 

0.3 

0.5 

372.8472 

y 

40 

0.1260 

0.1260 

0.3 

0.6 

- 

n 

- 

- 

- 

0.3 

0.7 

- 

n 

- 

- 

- 

0.3 

0.8 

- 

n 

- 

- 

- 

0.3 

0.9 

- 

n 

- 

- 

- 



r 2i 

Fval. 

Conv. 

Iter. 

Pi A 

P2 

0.4 

0.4 

- 

n 

- 

- 

- 

0.4 

0.5 

372.8472 

y 

40 

0.1260 

0.1260 

0.4 

0.6 

- 

n 

- 

- 

- 

0.4 

0.7 

- 

n 

- 

- 

- 

0.4 

0.8 

- 

n 

- 

- 

- 

0.4 

0.9 

- 

n 

- 

- 
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Template of digit 3 with 13 landmarks 



Original Euclidean x-coordinates of the template 


Figure 2: Configuration template. 



r 2 

' 1 i 

' 2i 

Fval. 

Conv. 

Iter. 

p 2 

p -? 2 

0.6 

0.6 

372.8472 

y 

34 

0.1259 

0.1260 

0.6 

0.7 

372.8472 

y 

37 

0.1260 

0.1260 

0.6 

0.8 

- 

n 

- 

- 

- 

0.6 

0.9 

- 

n 

- 

- 

- 


There is no doubt that the algorithm does not work and cannot reach any credible solution; in 


fact for some starting points the algorithm does not converge even at a very low truncation of 10th 
degree. 


Meanwhile the polynomial density gives a fast convergence no matters the starting point we 
consider, the result in this case by different classical methods for optimization, is: 'p 2 = 0.9542 and 
p 2 2 = 0.6740. 


Now, if we want to approximate the exact result by using the result of Muihead, the series must 


be truncated in a large number as the following table shows. We highlight that the algorithm was 
intervened many times in the starting point in order to obtain an stable solution. 


Truncation 


ps 


Truncation 

pi 

P* 


Truncation 

Pi 2 

p? 2 

20 

0.5890 

0.5889 


120 

0.9246 

0.6697 


220 

0.9484 

0.6733 

40 

0.7319 

0.7318 


140 

0.9330 

0.6710 


240 

0.9500 

0.6735 

60 

0.8554 

0.6713 


160 

0.9389 

0.6719 


260 

0.9512 

0.6736 

80 

0.8922 

0.6670 


180 

0.9431 

0.6725 


280 

0.9521 

0.6736 

100 

0.9121 

0.6682 


200 

0.9461 

0.6730 


300 

0.9527 

0.6737 

Truncation 

pi 

P2 


Truncation 

Pi 

P 2 2 




320 

0.9531 

0.6737 


420 

0.9538 

0.6737 




340 

0.9534 

0.6737 


440 

0.9538 

0.6737 




360 

0.9536 

0.6737 


460 

0.9539 

0.6736 




380 

0.9537 

0.6737 


480 

0.9539 

0.6737 




400 

0.9538 

0.6736 


500 

0.9539 

0.6737 





With the implemented routines based on infinite series of zonal polynomials was impossible to 
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reach the solution obtained with the exact distribution. After extremely large truncation of 500 
they could not stabilize to the solution obtained with exact polynomial of 10-th degree. Note, that 
we cannot provide an analytical of numerical relationship between the truncation and the estimates, 
it is exactly the open problem addressed by Koev and Edelman (2006). 

Once the computation problem is solved, we must pay attention to the interpretation of the 
squared canonical correlation estimates, pf 2 = 0.9542 and p 2 2 — 0.6740, corresponding to the 
comparison between the template and a typical handwritten sample of digit 3. The first approach 
consist of computing the probability that a given handwritten digit 3 with particular rfi and r 2 2 
is far from the template, in this case we just need to compute -P((ri,r 2 ) > (r^, r 2 2 )), which is a 
feasible task by using the polynomial densities of Theorem |3.1| 

Now, to obtain a measure for interpretation, we can apply the same procedure to a non-template 
object, for example, any handwritten digit 3 selected at random, see Figure [3] the next table 
summarizes the convergent problem of optimization when we use the infinite series. First we 
demonstrate the inefficient of the series which is truncated at the level of the degree 10-th of 
the exact polynomial distribution. 


A given handwritten digit 3 



Original Euclidean x—coordinates 


Figure 3: A randomly selected handwritten digit 3. 
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0.5 

0.5 

0.5 

0.5 

0.5 


0.7 

0.7 

0.7 

0.8 

0.8 

0.9 


Fval. 


0.5 

0.6 

0.7 

0.8 

0.9 


553.8395 


Conv. 


Iter. 


36 


Pi 


0.1282 


Fval. 


Conv. 


0.7 

0.8 

0.9 

0.8 

0.9 

0.9 


Iter. 


Pi 


P 2 


P2 


0.1282 



2 

r 2i 

Fval. 

Conv. 

Iter. 

P* 

Pi 

0.6 

0.6 

- 

n 

- 

- 

- 

0.6 

0.7 

553.8395 

y 

40 

0.1282 

0.1282 

0.6 

0.8 

- 

n 

- 

- 

- 

0.6 

0.9 

- 

n 

- 

- 

- 


Again the malfunction of the algorithms is revealed in the simple truncation 10; no convergence is 
obtained for every starting point, in particular for the largest ones which are near to the mean. Each 
pair of solutions is so far from the estimates based on the exact density; those values were derived 
by different numerical methods, no problems of convergence was occurred and were independent of 
the starting points, the corresponding estimates are: pi = 0.9753 and p 2 2 = 0.8748. 

If we want to pursuit that solution we require enormous values for truncating the series, as it is 
shown in the following table. 


Truncation 


20 

40 

60 

80 

100 


Truncation 


320 

340 

360 

380 

400 

420 

440 

460 


Pi 


P 2 


0.5976 0.5976 

0.7446 0.7446 

0.8107 0.8107 

0.8481 0.8480 

0.8719 0.8718 


Pi 


P 2 


0.9642 

0.9658 

0.9671 

0.9682 

0.9692 

0.9701 

0.9709 

0.9715 


0.8721 

0.8724 

0.8725 

0.8726 

0.8727 

0.8727 

0.8728 

0.8728 


Truncation 


Pi 


P 2 


120 

140 

160 

180 

200 


0.8883 0.8882 

0.9188 0.8728 

0.9315 0.8708 

0.9400 0.8704 

0.9462 0.8706 


Truncation 


Pi 


Pi 


480 

500 

520 

540 

560 

580 

600 


0.9721 

0.9726 

0.9731 

0.9735 

0.9738 

0.9741 

0.9743 


0.8729 

0.8729 

0.8729 

0.8729 

0.8729 

0.8729 

0.8729 


Truncation 


Pi 


220 

240 

260 

280 

300 


0.9510 

0.9548 

0.9578 

0.9603 

0.9625 


P'2 


0.8709 

0.8711 

0.8714 

0.8717 

0.8719 


After a truncation of zonal polynomials of order greater than 600 the algorithm fails and no 
computational chance of reaching the estimation based on the exact polynomial distribution, is 
possible. 

Summarizing we have arrive at two estimates: Template vs Handwritten digit 3, p\ = 0.9542 
and p 2 2 = 0.6740; meanwhile Given Handwritten digit 3 vs Handwritten digit 3, pi = 0.9753 and 
p 2 2 = 0.8748, the difference is notable, but a way for comparison can be the addressed distribution of 
Theorem|3.1| In this case we just compute the probability that an extreme correlation structure as 


r 2 = 0.9753 and r 2 = 0.8748 due to a particular given handwritten 3, occurs under the distribution 
of the co rrela tion structure influenced by the template. In this case, integration of distribution of 
Theorem 3.1 provides that -P((r 2 ,r 2 ) > (0.9753,0.8748)) < 0.0015, which means that the Eulerian 
shape of the given handwritten digit 3 is so far from the usual template of number 3. Figure [4] 
shows the associated distribution of the template and an approximate region for integration. 


4.6 An application in machine vision 

We end this section by studying the Eulerian shape differentiation between two populations, the 
application can be useful in machinery vision. Consider two samples of 17 crackers as shown in 
figures [5] (cracker type A) and [^cracker type B). 

In cracker type A, we register 9 landmarks by hand, the five evident marks on the cracker and 
the four corners, see Figure [7] We try to place similar landmarks on cracker type B, corresponding 
to the four corners, the center and the four nearest holes to the corners, , see Figure [8] It is clear 
the notorious variance of placing this landmarks on cracker B versus the more easily registration 
on cracker A. For a complete process in machine vision, we can also use automatic corner detection 
as Harris’ method for the landmark registration, but we want to increase intentional the difficulty 
of the pursuit population discrimination. Note also the adequacy of objects for implementing the 
Eulerian differentiation, the hgures have strong symmetry and some flexibility of choosing the 
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Figure 4: Eulerian polynomial joint distribution of template vs handwritten digit 3 and a 
test for a randomly selected handwritten code. 



Figure 5: A 17-sample of crackers type A randomly selected, original objects. 


landmark labels is allowed because the geometric invariants of such shape; meanwhile, in affine or 
similarity studies, the labeling have more radical constraint. 

We ask for a significant difference in Eulerian shape between crackers A and B. Next we list the 
results of some approaches. 
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Figure 6: A 17-sample of crackers type B randomly selected, original objects. 



Figure 7: A 17-sample of crackers type A randomly selected with 9 landmarks. 


METHOD 1: USING A DEGENERATED TEMPLATE: 

We proceed with a similar methodology of the handwritten digit 3 landmark data. In this case 
the selected template is just a rectangle with its center, in such way that 5 of the 9 landmarks 
are concentrated in the geometrical center and the remaining ones are the vertices, see Figure [9] 
Corresponding estimates for population of cracker A are pi = 0.7959 and p 2 2 = 0.7955; meanwhile 
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Figure 8: A 17-sample of crackers type B randomly selected with 9 landmarks. 


estimates for population of cracker B are, pi = 0.6218 and p 2 2 = 0.6212, as we expect, the template 
is a little more correlated with cracker A than cracker B. In this case we just compute the probability 
that an extreme correlation structure as r 2 = 0.7959 and rf = 0.7955 due to population correlation 
structure estimates of cracker B with the given template, occurs under the distribution of the 
correlation structure infl uenc ed by the template and population of cracker A. After integration of 

we have that -P((r 2 ,r 2 ) > (0.7959,0.7955)) < 0.0141, which means 


3.1 


distribution of Theorem 

that the Eulerian shape of cracker B is different from the Eulerian shape of cracker A, even under 
such degenerated template with more than 50 percent of coincident landmarks. Figure [10] shows 
the associated distribution of the template and an approximate region for integration. 

METHOD 2: USING A SIMILAR TEMPLATE TO CRACKER A: 

If we select a template where 4 of the above coincident landmarks are equally spaced in the 
diagonals of the rectangle (very similar to cracker A), see Figure 11 the corresponding estimates 
for population of cracker A are pf 2 = 0.9751 and p 2 2 = 0.9747; meanwhile estimates for population 
of cracker B are, f>i = 0.9714 and p 2 2 = 0.9710. Then we compute the probability that an extreme 
correlation structure as r 2 = 0.9751 and r\ = 0.9747 due to population correlation structure 
estimates of cracker A with the given template, occurs under the distribution of the correlation 
structure i nflue nced by the template and population of cracker B. After integration of distribution 
of Theorem 3.1 we have that P((rf, r 2 ) > (0.9751,0.9747)) < 0.0017, which means that the Eulerian 
shape of cracker A is so far from the Eulerian shape of cracker B. Figure [l2| shows the associated 
distribution of the template and an approximate region for integration. 

METHOD 3: USING A TEMPLATE CLOSE TO CRACKER B: 


If we select a template where 4 of the 5 coincident landmarks of method 1, are place on the sides 
of the rectangle (more closer to cracker B), see Figure[l3]the corresponding estimates for population 
of cracker A are pf 2 = 0.8928 and p 2 2 = 0.8927; meanwhile estimates for population of cracker A 
are, pi = 0.9569 and p 2 2 = 0.0.9567. Then we compute the probability that an extreme correlation 
structure as r 2 = 0.9569 and r 2 = 0.9567 due to population correlation structure estimates of cracker 
b with the given template, occurs under the distribution of the correlation structure influenced by 
the template and population of cracker A. After integration of distribution of Theorem |3.1| we have 
that -P((r 2 ,r 2 ) > (0.9569,0.9567)) < 0.00038, which means that the Eulerian shape of cracker A 
is so far from the Eulerian shape of cracker B. Figure |~i~T| shows the associated distribution of the 
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Figure 9: A degenerated template 



Figure 10: Eulerian polynomial joint distribution of a degenerated template vs cracker B 
and a test for a cracker A. 


template and an approximate region for integration. 

Note that the three methods ratifies the same conclusion about Eulerian shape difference, which 
seems to be independent on the template (of course it must satisfies the geometrical invariants 
demanding for Eulerian shape). In some experiments is difficult to suggest a template, but the 
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Figure 11: A template close to cracker A. 


tn 1»2j 2*2| and P<[l 1 1 2,r 2*2)>|0.975l.0.9747i)| 



1 1 


Figure 12: Eulerian polynomial joint distribution of a cracker B and template similar to 
cracker A, and a test for cracker A. 


addressed Heuristic independence shows that it is not a matter of big difference among the two 
pairs of estimates, but the joint distribution tends to be so narrow when the estimates are close as 
in method 2, meanwhile a degenerated template will provide apart enough the pairs of estimates 
but the joint distribution for testing is so wide, as we can see in method 3. Medium difference as 
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Figure 13: A template close to cracker B. 



Figure 14: Eulerian polynomial joint distribution of a cracker A and a template similar to 
cracker B, and a test for cracker B. 


in method 1 has a proportional joint distribution size. The most important fact is that in presence 
of shape differences, the selected template emerges as a theoretical or artificial object to perform 
the test, but it is irrelevant (in the sense of Eulerian restrictions). A very far artificial template will 
separate differences but the associated distribution has big variance, and it affects the resolution 
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of the test. Some mixing ideas for an artificial template close to intuition could be a good choice. 
Another usual conjecture for the template is a real observation of the sample, however, we have to 
be sure that the assigned template really belongs to the addressed population. 

We try to apply this methodology with a very sensitive landmark data, the schizophrenia land¬ 
mark date due to Bookstein (1996). In this case 13 landmarks in 2 dimensions, were register on 
a near midsaggital brain scan of 14 controls and 14 patients with schizophrenia. The landmarks 
choosen in the MRI were: (1) splenium, posteriormost point on corpus callosum; (2) genu, anteri- 
ormost point on corpus callosum; (3) top of corpus callosum, uppermost point on arch of callosum 
(all three to an approximate registration on the diameter of the callosum); (4) top of head, a point 
relaxed from a standard landmark along the apparent margin of the dura; (5) tentorium of cerebel¬ 
lum at dura; (6) top of cerebellum; (7) tip of fourth ventricle; (8) bottom of cerebellum; (9) top of 
pons, anterior margin; (10) bottom of pons, anterior margin; (11) optic chiasm; (12) frontal pole, 
extension of a line from landmark 1 through landmark 2 until it intersects the dura; (13) superior 
colliculus; see Bookstein (19961 and Dryden and Mardia ( 1998| ). 

First we use as a template one of the MRI sampled from schizophrenic patients, see Figure |15| 
The computations are simpler and no new integrals must be computed in Theorem |3.1[ because this 
data base has the same number of landmarks than the previously studied problem of handwritten 
digit 3. Again the distribution of the squared canonical correlation population is just a polynomial 
of degree 12 in the latent roots of the matrix argument. 


The corresponding estimates for population of schizophrenic patients are pi 2 = 0.9790 and p 2 2 = 
0.9789; meanwhile estimates for population of controls are a little less as we expect, pf = 0.9783 and 
P 2 = 0.9781. Then we compute the probability that an extreme correlation structure as rf = 0.9790 
and rf = 0.9789 due to population correlation structure estimates of schizophrenic patients with 
the given template, occurs under the distribution of the correlation structure influenced by the 
template and population of controls. After integration of distribution of Theorem |3.1| we have 
that P({ri,r 2 ) > (0.9790,0.9789)) < 0.0024, which means that the Eulerian shape of schizophrenic 
patients is different from the Eulerian shape of controls. Figure [16] shows the associated distribution 
of the template and an approximate region for integration. 



C 0 


Figure 15: A template corresponding to a randomly sampled schizophrenic patient. 


A similar result is obtained by taking an arbitrary template randomly selected from the sample 
of controls. 
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Figure 16: Eulerian polynomial joint distribution of controls and a schizophrenic patient- 
template, and a test for schizophrenic patients. 


Concluding remark: Finally, observe that the exact Eulerian polynomial joint distribution can 
be used for testing the structure of correlation, as in the same way that equality mean shape tests 
are usually done in shape analysis, with the difference that we are using the exact distribution 
instead of the asymptotic distributions of the statistics traditional implemented. Of course, the 
inference on correlation structure differs from the inference about mean shape, and both analysis 
complement each other, but it seems that the Eulerian shape analysis is sufficient for detecting 
differences between populations as the similar mean shape analysis. In fact, we have applied the 
methods presented here to some of the classical landmark data studied in literature by asymptotic 
distributions or inference in the tangent plane as those presented in Dryden and Mardia (19981, 
and we have obtained the same conclusions about mean shape differences, among the applications 
we do not include in this manuscript we count the following landmark data: mouse vertebra, gorilla 
skulls and macaque skulls. However, for completeness, a subsequent work will include the addressed 
associated Eulerian mean shape theory and we will back to this applications. 

Observe, that the exact Eulerian polynomial joint distribution can be used for testing the 
structure of correlation, as in the same way that equality mean shape tests are usually done in 
shape analysis, with the difference that we are using the exact distribution instead of the asymptotic 
distributions of the statistics traditional implemented. Of course, the inference on correlation 
structure differs from the inference about mean shape, and both analysis complement each other. 
In fact, a subsequent work will include the addressed associated Eulerian mean shape theory and 
we will back to this applications. 
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